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ABSTRACT 

We propose a method for measuring the cosmological density parameter from 
the statistics of the expansion scalar, 6 = if _1 V -v, - the divergence of peculiar velocity, 
expressed in units of the Hubble constant, H = lOOh km s _1 Mpc _1 . The velocity field is 
spatially smoothed over ~ 10h -1 Mpc to remove strongly nonlinear effects. Assuming 
weakly-nonlinear gravitational evolution from Gaussian initial fluctuations, and using 
second-order perturbative analysis, we show that (# 3 ) oc —Q~°' 6 (6 2 ^ . The constant 
of proportionality depends on the smoothing window. For a top-hat of radius R and 
volume-weighted smoothing, this constant is 26/7 — 7, where 7 = — dlog (# 2 ) /dlog R. If 
the power spectrum is a power law, P(k)ock n } then 7 = 3 + n. A Gaussian window yields 
similar results. The resulting method for measuring is independent of any assumed 
biasing relation between galaxies and mass. The method has been successfully tested with 
numerical simulations. A preliminary application to real data, provided by the POTENT 
recovery procedure from observed velocities favors ~ 1. However, because of an uncertain 
sampling error, this result should be treated as an assessment of the feasibility of our 
method rather than a definitive measurement of 0. 



keywords: Gravitation, Cosmology: theory - dark matter - large scale structure of the 
Universe 
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1. INTRODUCTION 

The standard model of large-scale structure assumes that the structure we 
observe today developed by gravitational instability from small-amplitude, random initial 
fluctuations, induced in the early universe, perhaps during an inflation phase (see Peebles 
1993, for a review). The most natural hypothesis (see Bardeen et al. 1986 and 
references therein) is that the initial probability distribution of the density contrast, 
8(x) = p(x)/(p^ — 1, was Gaussian (here p is the matter density, x are the comoving 
coordinates, and ( • • • ) denote ensemble averaging). Under this hypothesis, the initial 
distribution has zero skewness, (8 3 ^ = 0. With the growth of gravitational clustering, 
deviations from the initial Gaussian distribution develop. In particular, to second order in 
perturbation theory, the gravitationally induced skewness and variance are related by the 
scaling 

the above result assumes (S 2 ) < 1 and O = 1 (Peebles 1980, §18). Bouchet et al. 
(1992) derived S3 for arbitrary O. The O-dependence turned out to be extremely weak, 
in agreement with earlier suggestions based on numerical integration of the equations of 
motion (Martel & Freudling 1991). To mimic the observational techniques of estimating 
the moments of 6, Juszkiewicz et al. (1991;1993a) derived the S3 parameter for a spatially 
smoothed field. They generalized earlier results (Grinstein & Wise 1986, Goroff et al. 
1986), and showed that the filtering introduces a dependence of S3 on the slope of the 
power spectrum of the density fluctuations, while the weakness of the O dependence 
of S3 is preserved. These perturbative results were successfully tested against N-body 
simulations (Weinberg & Cole 1992; White 1992; Bouchet & Hernquist 1992; Juszkiewicz 
et al. 1993a, b; Lahav et al. 1993; Lucchin et al. 1994, Bernardeau 1993). 

The weak dependence on a such a poorly known parameter, as O, makes 
observational estimates of S3 useful in constraining the initial probability distribution of 8 
and its intrinsic (as opposed to the gravitationally induced) skewness (Silk & Juszkiewicz 
1991; Jaffe 1993, Bartlett et al. 1994). The observed density field of galaxies, under 
a reasonable assumption of local biasing relation between the galaxy density and the 
underlying mass density, was indeed found to be consistent with the scaling relation (1) 
(Coles & Frenk 1991; Bouchet et al. 1991; 1993 - for the IRAS survey; Gaztanaga 1994 
for the APM survey; for the effect of local biasing on S3, see Fry & Gaztanaga 1993 
and Juszkiewicz et al. 1993b). These results are consistent with the Gaussian initial 
conditions hypothesis and imply that the skewness present in the data seems to be entirely 
induced by gravity. Similar conclusions were reached by Nusser et al. (1993;94), who used 
the Zel'dovich (1970) approximation to reconstruct the initial probability distribution of 
8 from the IRAS data, and by Feldman et al. (1994), who investigated the spread of 
fluctuations in the measurements of the IRAS power spectrum. 

The independence of S3 on O, so useful for the above purposes, makes it useless 
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for measuring 0. However, as we will show in this paper, there is a close relative of S3, 
which can be used to measure the cosmic mass density parameter. This is the skewness of 
the distribution of the divergence of the peculiar velocity (see Dekel 1994 for a review on 
cosmic flows). Radial peculiar velocities of more than 3000 galaxies have been measured 
so far within a sphere of radius ~ 7000 kms -1 about the Local group, and the POTENT 
reconstruction method enables the processing of this data into a three-dimensional velocity 
field, smoothed at ~ 12 h _1 Mpc, based on the asserted potential-flow nature of quasilinear 
gravitating flows. This irrotational velocity field, v(x), can be characterized, up to an 
additive bulk velocity, by a single scalar field such as the velocity divergence 

6{x) = -^V-v(x) , (2) 
H 

called the expansion scalar (see, e.g. Peebles 1980). Division by the Hubble constant, H } 
makes 6 dimensionless. 

The gravitational instability equations governing the evolution of fluctuations of a 
pressureless gravitating fluid in a standard cosmological background during the matter era 
are the continuity equation, the Euler equation of motion, and the Poisson field equation 
{e.g. Peebles 1993): 

6 + V-v + V-(v6) = , 
v + 2Hv + (v-V)v = -V$ , (3) 
V 2 $ = (3/2)H 2 QS . 

Here x. v, and $ are the position, peculiar velocity and peculiar gravitational potential in 
comoving distance units, corresponding to ax, <ir. and a 2 $ in physical units, with a(t) the 
universal expansion factor, t is the cosmological time and dots represent time derivatives. 
To linear order in 6, v and $, the growing-mode solution for the density fluctuation has a 
universal time dependence, 

6 1 (x,t) = D(t)e(x) . (4) 
The linear velocity is irrotational, derived from a velocity potential, 

v 1 {x,t) = -D(t)V(p(x) , (5) 

where tp is time-independent and obeying the linearized Poisson equation V 2 (,c = e, i.e. the 
corresponding linear gravitational potential is $1 = (3/2)HQ,D<p. The relation between 
the linear solutions of velocity and density in the linear approximation is then 

e 1 (x) = -f(Q)S 1 (x) , (6) 

where /(O) = if -1 D/D O ' 6 (see Peebles 1993, eq. 5.120). This type of relation has been 
used in many different ways to determine O from the comparison of observed velocity and 
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density fields (see Dekel 1994 §8 for a review), but since the density field is derived from 
the distribution of galaxies, such determinations are necessarily subject to an arbitrary 
choice of biasing relation between the densities of galaxies and mass. 

In the present paper we derive from the velocity field alone: our method is 
independent of galaxy biasing. The linear relation (6) indicates that a relation between the 
moments of analogous to equation (1) should depend strongly on and may therefore 
be useful for measuring 0. Equation (6) also indicates that the perturbative technique 
used to compute the moments of 8 should be applicable in the well. This 

calculation is presented in §2 for a general value of 0, taking into account smoothing with 
a Gaussian or a top-hat window function and an assumed power spectrum. The results, 
derived by perturbative calculations with respect to the initial fluctuations, can only be 
applied to a heavily-smoothed velocity field, where the fluctuation amplitudes are not too 
nonlinear. The validity of the perturbative approximation for deriving under realistic 
nonlinear conditions is confirmed using numerical simulations in §3, where we also present 
preliminary results based on real observational data. Our method is discussed in §4. 

2. MOMENTS OF THE EXPANSION SCALAR 

2.1. Unsmoothed Fields 

We wish to compute the relation between the skewness and variance of 6, in analogy 
to equation (1), under gravitational evolution from Gaussian initial fluctuations. The 
perturbative expansion for is 

B = 1 +0 2 + O(8\) , (7) 

where 0\ is of order 8i, and 02 = 0(8\). The variance is then (0 2 ) = (0\) + 0(8^), so the 
leading term is simply 

<6» 2 >=/(fl) 2 <*?>. (8) 

The third moment, by equation (7), is 

(0*) = (01) + 3(0*02) + 0(8*) . (9) 

The skewness of 0\ vanishes as 0\ remains a Gaussian variable at all times, because it has 
a universal time dependence like 8\ (eqs. 4, 6). This leaves 3(0*02) as the leading term 
in equation (9), which implies that the perturbative expansion of must be carried out at 
least to second order. 

In calculating 02 we are guided by the analogous analysis of the density field (e.g. 
Bouchet et al. 1992; following Peebles 1980 §18), which yielded for the second-order term 



82(x,t) = D(ty[ 1 - 



F(t) 
1 + — — 
D(t) 2 



e 2 (x) + Ve(x)-V^(x) + ^ 



F(t) 

D(ty 



, a p(x)ip , a p(x) 

(10) 
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The subscript ", a/3" denotes spatial derivatives with respect to the a and f3 components 
of x. In addition to a growing mode oc D(t) 2 , the second-order solution has a growing 
mode ocF(t), where F(t) is the solution of the equation 



F + 2HF - -H 2 QF = -H 2 QD 2 . 
2 2 



(11) 



The ratio [F(t) / ' D(t) 2 \ turns out to be nearly independent of time, which implies that it 
is nearly independent of (Bouchet et al. 1992), 

|r - ?n- 2/M - (12) 

The corresponding second-order term of 6 is then obtained by substituting 8\ + 82 and V\ 
into the continuity equation and equating the second order terms, 

e 2 (x,t) = -/(0)D 2 {Ce 2 (x)+Ve(x)-V^(x) + [l-C]^, a ^)^, a ^(x)}, (13) 

F 



C(0) 



(14) 



2DD 

Much like F/D 2 in eq. (12), we find in a similar way that C is also almost independent of 
L>. 



c(o)^^- 1 / 21 . 



(15) 



This approximation is accurate to within 2% in the range 0.1 < O < 10 (Figure 1). A 
better, empirical approximation, accurate to within 0.5% in the range 0.1 < O < 1, is 
given by C(fi) w 3/7 O" - 055 . 



0.5 




Fig 1. The function C(O) of Eq. (14) (solid) compared to the approximate forms 3/7 ' 
(dashed), and 3/7 °- 055 (long dashed). 
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Substituting equation (13) in the second term of equation (9) we obtain (e.g. 
following Peebles 1980, §18, but for a general O) 

(e 3 ) = -m)[AC + 2] (o^) 2 ) 2 , (16) 

and then, with equation (8), 

This is a key result, valid at the lowest order in (# 2 ). As expected, the O dependence of 
T3 is indeed quite strong, which promises that T3 can provide a sensitive measure of O. 

2.2. The Effects of Smoothing 

The ratio in equation (17) refers to the statistics of the unsmoothed 6(x) at a point. 
In practice, however, the observed velocity field must involve significant smoothing of 
fluctuations at the small-scale end, in order to overcome inevitable noise and non-linear 
effects. The smoothing is expected to affect T3, as it affects S3. 

The velocity field smoothed by a window function W(x) is defined by a convolution 
in real space, 

v s (x) = J d 3 x'v(x') W(x - x') . (18) 

The divergence of the smoothed velocity, V-u„, equals the smoothed divergence of the 
velocity, H8 S} because the convolution and the divergence commute as linear operators. It 
can therefore be written as the Fourier integral 

6 s (x) = (2tt)- 3 / 2 J d 3 k W(k) 6 k e ik x , (19) 

where W and 6 are the Fourier transforms of W and 6 respectively. The variance, to 
leading order, is then 

(6 2 s ) = £> 2 (2tt)- 3 J d 3 kP(k)W(kf . (20) 

We used the first order relation 6(h) = De(k) } where e(k) is the Fourier transform of e(x), 
and the power spectrum is 

{e{k)e{k')) = 6 dirac (k + k')P(k) . (21) 
The skewness of the smoothed divergence becomes 

<^> = -f(Q)D 4 (2Tr)- 6 J d 3 Jb d 3 Jb' P{k) P(k') W{k) W(k') W(\k + k'\) T(ifc, k') , (22) 
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T(k, k') = 6C + Q^k/k') + 6 [1 - C] i? , (23) 

with \i = (k ■ k')/kk'. For a power-law spectrum, P(k) oc k n , an integral similar to the 
one above was reduced to an analytic expression by Juszkiewicz et al. (1993a) both for 
a top-hat and for a Gaussian window, and it is straightforward to rederive T3 using their 
methods. Recently the calculation for a Gaussian window function was extended beyond 
the integer values of n, considered by Juszkiewicz et al. (1993a), onto arbitrary n > —3: 
Lokas et al. (1994), expressed the integral (22) in terms of the hypergeometric function. 
The resulting skewness parameter is given by 



n- 1 



3F ( Z±*, Z±*; 3 -; + (4C-n-4)F ( ± 
1 2 ' 2 ' 2' 4/ v ; V 2 ' 2 ' 2' 4 



(24) 



The hypergeometric function F should not be confused with the perturbation growth rate 
F(t). The results for integer n are listed in Table 1; Figure 2 shows T3 as a function 
of and n. We notice that the dependence of T3 on the type of the window becomes 
stronger for higher values of n, where the small-scale structure is more important. The 
dependence remains strong for all spectral slopes. 

Table 1. Skewness vs. spectral slope for Gaussian and top-hat smoothing windows, can be 

read off Figure f , or using the approximations discussed in the text (eq. f 5). The terms in curly brackets 
have to be added in a future case of galaxy- weighted smoothing, where 6 is the linear biasing factor of the 
sampled tracers. 

n — /(0)T3, Gaussian window function — /(fi)T3, top-hat 

-3 2 + 4C(fi) 2 + 4C(fi) 

-2 6 (f - f ) + 6C(fl) (V3 - f ) - {66 (l - f ) } 1 + 4C(fi) - {b} 

-1 6 - 6) + 6C7(fl) (8 - J|) - {66 - l) } 4C(0) - {26} 

6 (f - ^) + 6C(fl) - f ) - {66 (l - ^) } -1 + 4C(fi) - {36} 

1 6 (i% - !) + (! - - { 6& ~ } - 2 + " W 



Bernardeau (1993) showed that for the top-hat window, the 'pure power-law' results 
of Juszkiewicz et al. (1993a) remain valid when n is not an integer and is varying with the 
smoothing radius, R: 

T - 1 ^4-9 ! - dlQg(^) 

T3 = ^y = -m (iC+2 lh i= -^f- (25 > 
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Here the spectral index is allowed to vary in the range > 7 < 4; outside of this range 
T3 diverges. In the top-hat case, 6 S represents a local fluctuation in the mean expansion 
(because it is the integral of V-t> in a sphere of radius R } which, by the Gauss theorem, is 
proportional to the total flux of matter across its spherical boundary and thus represents 
a local perturbation in the Hubble flow). For a pure power-law spectrum the variance of 8 
in a top-hat sphere of radius R is (8 2 ^ oc R~( 3+n \ so 7 = 3 + n, and for a top-hat window 
we have 

Ts ~ (y ~ ( n + 3 )) ' valid for - 3 < n < 1 , (26) 

to a good approximation, based on ignoring the O-dependence of (7(0) and setting C ~ 3/7. 

Since the Zel'dovich (1970) approximation is frequently used to study 
weakly-nonlinear dynamics and statistics, it is instructive to rederive T3 from a 
perturbative expansion of the Zel'dovich solution. We will use the superscript 'Z', to 
distinguish the result from T3. The Zel'dovich solution gives 82, similar to that in eq. (10), 
but with F = (Grinstein & Wise 1987, Bouchet et al. 1992). Therefore, C = 0, and the 
term 26/7 in equation (26) is replaced by 2, giving 

Tf = -fi-°- 6 (n + 1) ^ T 3 . (27) 

Numerical simulations, discussed in the next section as well as those of Juszkiewicz et 
al. (1993b) agree with the real perturbation theory (T3) but not with Tf. Note, that 
for spectral slopes —1 < n < 5/7, the Zel'dovich approximation fails to reproduce even 
the sign of T3. This failure is not surprising since the perturbative solutions of the real 
equations of motion conserve mass as well as momentum to all orders by definition, while 
the Zel'dovich solution conserves momentum only to first order. 

We have also computed the correction for galaxy-weighted smoothing, which may 
have a potentially interesting future application. In this case the distribution of 6 is 
heavily weighted by densely-sampled regions. Since the sampling density of velocity tracers 
is correlated with the underlying galaxy density, which is probably correlated with the 
density of matter ("biasing") and thus with 6, the distribution of 6 smoothed is severely 
biased. In the ideal case where the density fluctuation of sampled galaxies is proportional 
to the mass-density fluctuation, 8 saTnp i e = b saTnp i e 8 } we find the smoothing correction to 

consist of an additional term in eq. (22), +66/ -1 (#«) 2 , and an additional term in eq. (23), 
+66(1 + k ■ k' /k' 2 ). For a top-hat window, the final expression replacing eq. (25) is 

T 3 = -j^[AC + 2-(l + b sample ) 1 ] . (28) 

Future data may have uniform, complete sampling, with b saTnp i e the same as 6, the linear 
biasing factor of all galaxies. Since the galaxy- weighted T3 depends on O ' 6 and on 
f3 = 6/0° ' 6 in different ways, it may be used to determine O and b separately. Methods 
based on comparing the velocity data to the galaxy distribution are limited to determining 
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the degenerate combination /3, so the galaxy- weighted T3 may eventually become an 
interesting tool. In Table 1, we list in curly brackets the extra terms, which should appear 
when galaxy- weighted smoothing is applied; the b parameter there denotes b saTnp i e . For 
a Gaussian window, the expression in square brackets in eq. (24) should contain the 
following additional terms: +6{3F(^± 1 , ^±1; §; f)-nF(^p, ^ f ; i)~6} • The results 
obtained in this way are significantly different from those we have just obtained for a 
volume-weighted smoothing. 




Fig 2. The exact second-order values of — T3 = — \@r)/\@r) as a function of O for different 
values of the spectral index Tl of a power-law spectrum (n = — 3, — 2, — 1,0 etc. from top to bottom). 
Left: smoothing with a Gaussian window. Right: smoothing with a top-hat window. The filled squares 
correspond to CDM numerical simulations with fi = 0.3 and fi = 1.0, Gaussian smoothed with radius 
12 h _1 Mpc (where n 0). The error bars represent the cosmic scatter due the finite-size sample of 
radius 40 h "'Mpc. The dashed lines correspond to approximation (26). 



3. PRELIMINARY APPLICATION TO SIMULATIONS AND REAL DATA 



3.1. Simulations 
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The validity of the second-order approximation for T3 in the case = 1 was 
confirmed using N-body simulations (Juszkiewicz et al. 1993b). In this paper we use similar 
simulations for two purposes: to test the dependence of the second-order approximation 
derived for T3, and to estimate the cosmic scatter due to the fact that, in practice, T3 is 
computed from one finite realization of a random field. This is only a preliminary test 
using simulations; the results are expected to be model dependent, sensitive, in particular, 
to the shape of the fluctuation power spectrum. 

Twelve simulations were run (Nusser & Dekel 1993) for each of the cases, = 1 
and 0.3, using a particle-mesh code (Bertschinger and Gelb 1991) with 64 3 grid cells and 
particles in a periodic cubic box of comoving size 160 h _1 Mpc. The initial conditions for 
each simulation were a random Gaussian realization of the "standard" CDM spectrum 
(Davis et al. 1985, h= 0.75). The simulations were stopped at a time step when, based 
on linear growth, the rms density fluctuation in top-hat spheres of radius 8 h _1 Mpc 
was unity (erg = 1). The velocity field was evaluated in a 64 3 cubic grid of spacing 
2.5 h _1 Mpc by linear cloud-in-cell assignment of velocity to a grid point from the particles 
in its neighboring cells. The fields were then smoothed on a much larger scale, volume 
weighted, using a spherical Gaussian window of radius 12 h _1 Mpc. This procedure ensures 
smoothing over nonlinear and orbit-mixing regions and it produces quasilinear fields which 
are comparable in rms amplitude to the fields deduced from observations. 

We select from the simulations many independent spheres of radius R = 40 h _1 Mpc, 
compute T3 within each sphere by summing over the grid points, and obtain the mean and 
standard deviation of T3 over the spheres. The mean is found to be T3 = —1.8, —4.1 
for = 1, 0.3 respectively, in pleasant agreement with our second-order predictions (Fig. 
2). The cosmic standard deviation at R = 40 h _1 Mpc, for the assumed power-spectrum, 
is ei = ±0.7, ±1.3 respectively. It becomes smaller at larger volumes, of course, but 
R ~ 40 h _1 Mpc is the practical current limit for high-quality data. 

3.2. Data 

In the present method-oriented paper, we bring only a preliminary application to 
real data, and defer a careful, detailed application to a paper in preparation. The 6 field 
is provided by the current POTENT analysis (Dekel et al. 1994; based on Mark III data, 
Willick et al. 1994). The raw data consist of ~ 3000 galaxies with measured redshifts 
Zi and redshift-independent estimated distances r{ of assumed random Gaussian errors 
G{ ~ (0.15 — 0.21)ri. The peculiar velocities are U{ = Z{ — r{. The distances are corrected 
for inhomogeneous Malmquist bias by heavy grouping (Willick 1991) and by adopting 
a correction procedure which assumes that IRAS 1.2 Jy galaxies trace the underlying 
galaxy distribution. The radial peculiar velocities are first smoothed with a Gaussian 
window of radius 12 h _1 Mpc onto a spherical grid, yielding a radial velocity field u(x). 
The smoothing suppresses the random uncertainties due to distance errors and sparse, 
nonuniform sampling. The adopted weighting scheme attempts to mimic equal-volume 
weighting and thus to minimize the associated systematic biases (Dekel, Bertschinger & 
Faber 1990, DBF). POTENT then recovers the transverse components of the velocity 
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field by requiring that it is derived from a scalar potential, v = — Vo, as predicted by 
linear gravitational instability and maintained later via Kelvin's circulation theorem under 
appropriate smoothing (Bertschinger & Dekel 1989). The potential at each point x is 
computed by integrating <f>(x) = — J Q r u(r', i9, cp)dr' , and the transverse velocity components 
are derived by differentiating this potential. The divergence field, 6, is computed by finite 
differencing using a cubic grid of 5 h _1 Mpc spacing. 

The resulting divergence field is contaminated by random errors, 62, due to the 
uncertainty in the distance measurements, 0{. They are estimated at each grid point by 
Monte Carlo noise simulations (DBF). The errors at r = 40 h _1 Mpc are typically smaller 
than €2 = 0.2. These random errors should not have a systematic effect on the mean 
and the skewness of 6, but they do add systematically to the true variance. We therefore 
compute the variance of the Monte-Carlo errors within the 40 h _1 Mpc sphere and subtract 
it from the total variance computed for the field recovered by POTENT. The same set of 
Monte-Carlo noise simulations allow us to estimate the distance error in T3. 

The value of T3 within 40 h _1 Mpc is tentatively found to be T3 = —1.1 ± 0.7. For 
the total error, the quoted distance error, €2 = ±0.7, should roughly add in quadrature to 
the cosmic scatter estimated by the simulations, e\. Hence, fi = 1 is fully consistent with 
the data, while fi = 0.3 is marginally rejected at the ~2-sigma level. 

4. DISCUSSION 

We sketched a method for measuring fi from the probability distribution function 
(PDF) of the smoothed velocity divergence field under the assumption of mildly-nonlinear 
gravitational growth from Gaussian initial fluctuations. Unlike the PDF of 6, whose 
shape in the quasi-linear regime is uniquely determined by its variance, the PDF of 6 
at a given variance depends on the dynamics which made it grow, via fi. In particular, 
the ratio of moments T3 = (# 3 )/(# 2 ) is expected to be oc fi -0 ' 6 . The expected value 
of T3 was computed using second-order perturbative approximation and confirmed by 
N-body simulations. The result takes into account the inevitable effect of smoothing, 
which depends on an assumed power spectrum. 

When applied to the volume- weighted velocity field recovered by POTENT from the 
observed radial peculiar velocities of galaxies, the method is free of assumptions regarding 
the biasing relation between galaxy density and mass density. The sampled galaxies are 
assumed to trace the large-scale velocity field, like any other test body subject to the 
large-scale gravitational field and sharing the same initial conditions, but they are not 
assumed to trace the underlying mass distribution. This method is thus ascribed to a 
unique family of methods which give constraints on fi directly (including, e.g. , Nusser & 
Dekel 1993; Dekel & Rees 1994), unlike most methods which measure the degenerate ratio 
fi°- 6 /& (see Dekel 1994 for a review). 

The method relies on the hypothesis that the initial fluctuations were a realization 
of a Gaussian random field. As mentioned in the Introduction, this hypothesis finds 
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support in the observed IRAS galaxy density, but the hypothesis can in fact be tested 
using the velocity data itself. The scaling relation between the skewness and the square 
of the variance, (# 3 ) oc (0 2 ^ , is unique to Gaussian initial conditions. This relation can in 
principle be measured in the data by varying the smoothing length in the velocity analysis. 
A practical limitation arises from the limited extent of the available data compared to the 
sparse sampling, which permits only a limited dynamical range for the smoothing length. 

In view of the data available, the method proposed is subject to two kinds of severe 
uncertainties. First, the random errors and biases in the recovered velocity field, due to 
the distance errors and the sparse nonuniform sampling, which we partly estimate using 
the POTENT Monte-Carlo noise simulations. Second, the fact that the data consists of 
only a few tens of independent points and therefore cannot be assumed a "fair sample" 
and cannot trace the PDF in great detail. The associated cosmic scatter is estimated using 
N-body realizations of the CDM spectrum, and similar estimates for other possible spectra 
would be useful. The two errors added in quadrature indicate that the method as applied 
to the existing data is capable of rejecting an = 0.3 cosmology with marginal confidence. 
With more data in a more extended volume the discriminatory power of the method will 
increase. 

The O dependence of the PDF of 6 is not limited to the skewness, and T3 is not 
necessarily the most useful measure of O. For example, Bernardeau (1993) computed 
a similar O-dependent function of moments involving the kurtosis of 6 (for top-hat 
smoothing). In practice, however, the skewness and kurtosis are both sensitive to the 
poorly-sampled, noisy tails of the distribution. Methods which make use of the whole 
shape of the PDF, that can be derived analytically in the quasi-linear regime (Bernardeau 
1994), might turn out to be more practical. 

Our preliminary results favor a high value of O, of the order of unity. This result is in 
accord with the result of Nusser & Dekel (1993), also based on the assumption of Gaussian 
initial fluctuations, combined with an analysis of the evolution of the PDF different from 
the analysis adopted here, and a different approximation for gravity in the quasi-linear 
regime - the Zel'dovich approximation. High O values are also obtained from the velocity 
divergence near voids (Dekel & Rees 1994), and by comparing the velocity data to the 
galaxy density (Dekel 1994 for a review). On the other hand, other estimates of O, such 
as those based on the anisotropy of galaxy clustering in redshift space {e.g. Cole et al. 
1994, and references therein) or variational principle analysis of trajectories of groups of 
galaxies (Peebles 1993), yield low O. 

It should be emphasized, that our results obtained from the data provide an 
illustration of the analysis rather than a definitive measurement of O because the present 
sets of data are too small. Our control over the sampling errors should improve in the next 
several years, when deeper surveys become available. 
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